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The BSSN and standard gauge equations are written in covariant form with respect to spatial 
coordinate transformations. The BSSN variables are defined as tensors with no density weights. 
This allows us to evolve a given set of initial data using two different coordinate systems and to 
relate the results using the familiar tensor transformation rules. Two variants of the covariant 
equations are considered. These differ from one another in the way that the determinant of the 
conformal metric is evolved. 



I. INTRODUCTION 

The BSSN formulation of Einstein's equations [TJ |5] is 
in widespread use in the numerical relativity community. 
These equations are most often used in conjunction with 
the "standard gauge" conditions, namely, 1+log slicing 
and the Gamma-driver shift. 

The BSSN variables include the conformal metric g a b, 
conformal factor and the trace and trace-free parts of 
the extrinsic curvature, K and A ab . They are defined in 
terms of the physical spatial metric g a b and the physical 
extrinsic curvature K a b by 

g a b = e 4v gab , (la) 
K ab = e^(A ab + g ab K/3) . (lb) 

The conformal metric is chosen to have unit determinant, 
g = 1, and A a b is trace-free. A key ingredient of the 
BSSN formulation is the use of the "conformal connection 
functions", defined by T a = —dbg ab - 

The 1+log slicing condition is an evolution equation 
for the lapse function a that takes the form [3] 

d t a = (3 c d c a - 2aK , (2) 

where f3 a is the shift vector. The Gamma driver shift is 
defined by @] 

d t /3 a = P c d c f3 a + ^B a , (3a) 

d t B a = p c d c B a + d t T a - f3 c d c T a - v B a , (3b) 

where B a is an auxiliary variable and rj is a constant. 
The term dtT a in Eq. (|3|d) is replaced with the right- 
hand side of the equation of motion for T a . 

The Gamma driver shift condition is not generally co- 
variant. In other words, Eqs. (J3j) do not preserve their 
form under a time-independent transformation of the 
spatial coordinates. To see this, note that (3 a and dt(3 a 
transform as contravariant vectors. The advection term 
(3 c d c f3 a in Eq. (|3^) does not transform as a contravariant 
vector. This spoils the covariance of Eq. ([3^). Note, how- 
ever, that the advection terms are not always included in 
the Gamma driver shift condition [5]. So for the mo- 
ment let us ignore the terms with (3 c d c acting on f3 a , B a , 
and T a . Then Eq. ^f) shows that B a and d t B a should 
transform as contravariant vectors. However, the right- 
hand side of Eq. d3p) depends on the connection functions 



r a , which do not form a contravariant vector. Instead, 
r a obeys a rather complicated transformation rule deter- 
mined from the following considerations. The conformal 
metric has unit determinant, g = 1. This equation is gen- 
erally covariant under spatial coordinate transformations 
only if g ab is defined as a type (°) tensor density of weight 
—2/3. This makes g a scalar, which is set equal to the 
scalar 1. Therefore the conformal connection functions 
r a = —dbg ab transform as the contraction of the deriva- 
tive of a type (p) tensor density of weight 2/3. This is a 
complicated transformation rule which I will not bother 
to write out in detail. The time derivative, dtT a , also sat- 
isfies this rule. We see that even if we ignore the advec- 
tion terms, there is a mismatch in Eqs. (J3j) in the way that 
the individual terms transform under time-independent 
changes of spatial coordinates. 

In addition to the conformal metric g a b and the confor- 
mal connection functions T a , the BSSN variables include 
ip, A a b, and K, Because g a b carries a nonzero density 
weight, the conformal factor ip must transform as the log- 
arithm of a weight 1 /6 scalar density. The variable A a b 
is a trace-free type (°) tensor density of weight —2/3. 
The trace of the extrinsic curvature K is a scalar. With 
these transformation rules, the BSSN equations and the 
1+log slicing conditions are covariant. 

In section II, I discuss the issue of general covariance 
in more detail. This serves as further motivation for the 
subsequent analysis. In Section III, I rewrite the BSSN 
equations in terms of simple tensor variables with no den- 
sity weights. One of the key steps in the analysis is the 
introduction of a background connection, as suggested 
by Garfinkle, Gundlach, and Hilditch [6 . Another key 
step is to recognize that the condition g — 1 should be 
replaced by an evolution equation for g [TJ. 

There are two natural choices for the evolution of 
g, which are presented in Section IV. One is the "La- 
grangian case" in which dtg = 0. Then g is constant 
in time, equal to its initial value. If the initial value 
of g is unity and the background connection vanishes, 
these equations reduce to the traditional BSSN equa- 
tions. These equations are not strongly hyperbolic unless 
the trace-free property of the variable A ao is actively en- 
forced during the evolution. As an alternative, one can 
add a term proportional to A to the evolution equation 
for the conformal metric. This yields a strongly hyper- 
bolic system without the need to actively enforce the con- 
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dition A = 0. 

Another choice for the evolution of g is d t g = Cpg, 
where C is the Lie derivative. This is the "Eulerian 
case". The BSSN equations in Eulerian form contain 
fewer terms than the traditional BSSN equations. 

In Section III the tensorial BSSN variables are used 
to rewrite the Gamma-driver shift condition in generally 
covariant form. If the initial value of g is unity and the 
background connection vanishes, the covariant Gamma- 
driver shift equations for the Lagrangian case reduce to 
the familiar Eqs. As an alternative, the shift con- 
dition can be defined using Eqs. ([3| with the Eulerian 
evolution equation for d t T a . 



II. SPATIAL COVARIANCE IN NUMERICAL 
RELATIVITY 

Let me discuss the issue of spatial covariance in con- 
crete terms. Let's say we are given a physical metric and 
extrinsic curvature, g ab and K abl that satisfy the Hamil- 
tonian and momentum constraints. For simplicity, let us 
assume that these tensors g ab and K ab are expressed in 
terms of a single coordinate patch with "Cartesian" coor- 
dinates x, y, and z. By calling the coordinates Cartesian 
I mean that each coordinate ranges over an interval of 
the real number line, with no periodic identification and 
no coordinate singularities. In this case we can trans- 
form to "spherical" coordinates r, 9, and <f> using the 
familiar relations x = rsin#cos</>, y — r sin (9 sin</), and 
z = rcos(9. Let us denote the initial data in spherical 
coordinates by g' ab , K' ab . The initial data in these two 
coordinate systems are related by 



9ab 

KL 



dx c dx d 

dx^dx^ 1 
dx c dx d 
dx^dx^ 



9cd i 

k cd 



(4a) 
(4b) 



where x a are Cartesian coordinates and x' a are spherical 
coordinates. 

We want to evolve these data using the BSSN system. 
Starting with the physical metric and extrinsic curvature 
in Cartesian coordinates, g ab and K ab , we apply the defi- 
nitions to obtain initial values for the BSSN variables 
g a b, P, A ab , K, and T a . Alternatively, we can start with 
the physical metric and extrinsic curvature in spherical 
coordinates, g' ab and K' ab , and define the BSSN variables 
g' ab , ip' , A' ab , K 1 , and T' a . The primed and unprimed 
BSSN variables will be related by the coordinate trans- 
formation rules outlined in the introduction. For exam- 
ple, the relation for the conformal metric is 



9ab 



dx 

dx* 



~ 2 ^ 3 dx c dx d 
dx^dx* 



9cd 



(5) 



The factor \dx/dx'\ is the Jacobian of the coordinate 
transformation. 



In order to evolve these data we must choose a lapse 
function a and shift vector (3 a . Let us assume for defi- 
nitencss that the lapse and shift are determined by evo- 
lution type equations and, for the moment, let us as- 
sume that these equations are spatially covariant. For 
example, one might consider 1+log slicing and a "modi- 
fied" Gamma driver condition obtained by replacing the 
term d t T a — f3 c d c T a with, say, D b Y. ab , and replacing the 
remaining spatial derivatives with covariant derivatives. 
(Here, D b is the physical spatial covariant derivative and 
S ab = (D a f3 b - aK ab ) TF is the distortion tensor. TF 
denotes the trace-free part.) 

These gauge conditions require us to specify initial val- 
ues for the lapse, shift, and auxiliary variable. Let a, 
(3 a , and B a denote these initial values in Cartesian co- 
ordinates. The initial values in spherical coordinates are 
related by the familiar tensor transformation rules: For 
the scalar lapse we have a! = a, and for the contravariant 
vector shift, 



dx' 



dx b 



73" 



(6) 



The auxiliary variable B a also transforms as a contravari- 
ant vector. 

Now we're ready to evolve both the "unprimed" Carte- 
sian coordinate data and the "primed" spherical coordi- 
nate data from the initial time ti to some final time tf. 
How do the primed and unprimed BSSN variables com- 
pare at tf? If, as we assumed above, the lapse and shift 
are determined by covariant relations, then the two sets 
of BSSN variables at time t f will be related by the same 
transformation rules that apply to the initial data. In 
particular, the two conformal metrics at time tf will be 
related by Eq. 

This, of course, is a good situation. We would like 
to have the option of evolving our initial data using dif- 
ferent spatial coordinate systems, and we would like to 
be able to compare the results. But the Gamma driver 
shift condition Eq. ^ is not covariant. Thus, at times 
t > t i} the shift vector obtained from the Cartesian co- 
ordinate evolution will differ geometrically from the shift 
vector obtained from the spherical coordinate evolution. 
These two shift vectors will not be related by a coordi- 
nate transformation; rather, they will be geometrically 
distinct vector fields. 1 

Because the shift vector depends on the coordinate sys- 
tem, the BSSN variables at time tf will not be related by 
the coordinate transformation rules outlined above. In 
particular, the conformal metrics at time tf will not be 
related by Eq. d5J. We can recombine the BSSN vari- 
ables to form the physical metric and physical extrinsic 



1 With the 1+log condition (jSl, the slicing of spacetime does not 
depend on the shift vector or the coordinate system. If the advec- 
tion term f3 c d c a is dropped, then a non— covariant shift condition 
will cause the slicing to depend on the coordinate system. 
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curvature. The physical tensors will not be related by 
the transformation Q. 

This shortcoming of the BSSN formulation with 
Gamma-driver shift can be fixed. This is done by rewrit- 
ing the BSSN equations and the standard gauge condi- 
tions in terms of regular tensors with no density weights. 
The density weights are removed by eliminating the re- 
quirement that the determinant of the conformal met- 
ric, g, should equal 1 [7j. It is then necessary to spec- 
ify an evolution equation for g. The simplest choice is 
to let g be constant in time. Then g is determined by 
its initial value, a weight 2 scalar density called g. The 
conformal metric g a b is then a type (°) tensor with no 
density weight. The next step is to define a variable 
A a = g bc (T bc — r£ c ) to take the place of the conformal 
connection functions [BJ. Here, T^ c are the Christoffel 
symbols constructed from the conformal metric, and T bc 
is a background connection. 

The BSSN and standard gauge equations can be writ- 
ten in terms of the regular tensors gab, ip, A a b, K, and 
A a . These equations reduce to the traditional forms in 
use by numerical relativity groups when the density g is 
set to unity and the background connection is that of a 
flat metric in Cartesian coordinates: T bc = 0. If we want 
to transform an initial data set from Cartesian to spher- 
ical coordinates, and preserve the tensor transformation 
rules under evolution, then we must transform g as a 
weight 2 density and T^ c as a connection. In particular, 
in spherical coordinates, g would no longer be 1 and T bc 
would no longer vanish. 

III. TENSOR VARIABLES FOR BSSN 

In this section I derive the BSSN equations from 
scratch using only tensors with no density weights. Be- 
gin with the "gdot-Kdot" form of the Einstein evolution 
equations: 

d±g a b = -2aK ab , (7a) 
d^Kab = akkab-2ak ac k c b 

+aR a b - D a D b a . (7b) 



The time derivative operator is defined by d± = dt — 

Cf3 where Cp is the Lie derivative along the shift. The 
Hamiltonian and momentum constraints are 

H = K 2 - k ab k ab + R , (8a) 

M a = D b (K b a - KS b a ) . (8b) 



Indices on the momentum constraint and extrinsic cur- 
vature are raised and lowered with the physical metric. 

The BSSN variables are defined by Eqs ([IJ. However, 
we will not assume any restrictions on the determinant 
of g a b or the trace of A a b- Then these definitions can be 
inverted to obtain 



9ab = (g/g) 1/3 g ab , (9a) 

Aab = (g/g)- 1/3 (Kab- g ab k/3 + g ab A/3) , (9b) 

V = 12 ln ^- 9 ) ' ( 9 °) 

K = K - A . (9d) 



Note that A = g ab A ab and K = g ab k ab - Defined in this 
way, the BSSN variables g a b, A a b, V, and K are tensors 
with no density weights. 

Now compute the time derivatives of Eqs. J9]) using 
the gdot-Kdot equations 0, then use Eqs. ([l]f to ex- 
press the results in terms of BSSN variables. This is a 
straightforward, although somewhat tedious calculation. 
It is useful to note that (g/g) 1 / 3 = e 4ip , and for the last 
two terms in Eq. ^p), 

Rab = Rab - 2D a D b ip + 4:D a (pD b ip 

-2g ab (D 2 ip + 2D c ipD c <p) , (10a) 
D a D b a = D a D b a ~ 4:D( a aD b )(p 

+2g ab D c aD c cp . (10b) 

The result of this calculation is 



d±g a b = - g a b d±ln g - 2aA a b+ 7^ag ab A , (11a) 
dj_A ab = ^A ab d_ L ]ng+~g ab d ± A-2aA ac Al + aA ab K + ~aA(5A ab - Ag ab - Kg ab ) 

+e-^[~2aDaDbV + 4aDavDbV + 4D (a aD b) ^-DaD b a + aRab] TF , (lib) 

d±<p = ~d±big- l -a{K + A) , (11c) 

d±K = -d±A + a(K + A) 2 + er iip (aR - &aD a tpD a (p - 8aD 2 ip - D 2 a - 2D a aD a f) . (lid) 

I 

The superscript "TF" denotes the trace-free part of the expression in brackets. 
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Now define 



AAL be — 1 be 1 be > 
AT a = g be AT a c ^ 



(12a) 
(12b) 



where T% c is a background connection. Although it is not 
necessary, it is convenient to think of the background con- 
nection as being built from a background metric, g ab . I 
assume that the background connection is time indepen- 
dent. Note that Ar a is a contravariant vector. 

From the definition of the Riemann tensor, we have 
the following identity: 

R ab = --g cd D c D d g ab + g c(a D b) AT c 



-9 Cd 9e(aRb)ed e + ff* AI^ AT {ab) 

+g cd (2ATl (a AT 



b)ed 



Ar^Ar 



<:bd 



(13) 



Here, D a is the covariant derivative and R a bcd is the Rie- 
mann tensor built from T% c . 



The equations of motion (111 imply 
<9j_(Ar a ) 



g bc D b D c p a - g bc R a bc<i (3 a 



2A ab d b a 



~D b (^A ab ) - 1 -AT a (d ± lng + 2aA) 

Vd 3 



g ab d b (d ± \ng-4aA) 



(14) 



Again, this is a straightforward, but tedious calculation. 
We now let A a denote a new BSSN variable which equals 
Ar a when the following constraint holds: 



C a = A a - AT a . 



(15) 



A a is a contravariant vector. 

Next, we modify the equations of motion using the 
Hamiltonian and momentum constraints. The Hamilto- 
nian constraint is TL — 0; from Eq. (|8^) we find 



n 



2 - { k + aY + \a 



2 1 1 a2 A ab A ab 



z- iv {R - 8D a <pD a ip - 8D 2 ip) . (16) 



Now add -aU to the right-hand side of Eq. JITJi). The 
momentum constraint is A4 a = 0; from Eq. d8p) we find 

M a = ^-e- iip b b {^fgA ab ) + 6e~ 4lp (A ab - Ag ab /3)d b <p 



V9 



-^g ab d b (2K/3 +A)+ e~ iv A bc AY a bc . (17) 



The equation of motion for A a is obtained by replacing 
dj_AT a on the left-hand side of Eq. (filf with dj_A a . We 



then add 2ae 4v M a to the right-hand side of this equa- 
tion. 

The analysis described above yields the following equa- 
tions for the tensor BSSN variables: 



dj_g a b 
d±A ab 

d±ip 
d ± K 

d±A a 



1 2 

^g a b d± lug - 2aA ab + -ag ab A , 

^A ab d± lug + ^g ab d±A - 2aA ac A c b + aA ab K + ^aA(5A ab - Ag ab - Kg ab ) 



+e -4v [-2aD a D b <p + 4aD a ipD b ip + 4D [a aD b) <p - D a D b a + aR ab ] 

— d±\ng- -a(K + A) , 
12 y 6 V ; ' 

-d x A + a (K 2 + 2KA) + aA ab A ab - (D 2 a + 2D a aD a ip) , 
3 

g bc D b D c (3 a - g bc R a bcd (3 d - ]-AT a 8 ± lug - ];g ab d b d ± lug 

3 6 

-2(A bc - g bc A/3)(6£d c a - 6aS^d c ip - aAT bc ) - \ag ab d b K , 



TF 



(18a) 

(18b) 
(18c) 
(18d) 

(18e) 



where 

T^-ab 



l -g cd DeD d 9ab + 9e(aD b) K c - g cd g e{a R b)cd e + g de AT c de AT {ab)c + g cd (2ATl (a AT b)ed + AY e ac AT ebd ^ 



(19) 



In Eq. (19 1, 7Z ab is defined by using A a in place of Ar a 



in the identity ( 13 1 



The Eqs. ( 18 ) are not complete evolution equations 



because d±g and d±A appear on the right-hand sides. 
These equations are consistent in the sense that we can 



use them to compute d±g, and the result is an identity: 
d±g = d±g. Similarly, Eqs. (18 1 yield an identity for 
d±A. 

Because the quantities d±g and d±A appear on the 
right-hand sides of Eqs. ( 18 ), we must specify how g and 
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A evolve. There are two natural choices for g, namely, 
d t g = and d±g = 0. In Ref. [7J [8], these were referred 
to as the Lagrangian case and the Eulerian case, respec- 
tively. For A, I will only consider the evolution equation 
dtA = 0. These cases are described in detail in the next 
section. 

Using the tensorial variables define above, the standard 
gauge conditions are: 



d t a = /3 c D c a - 2aK , 
d t /3 a = P c b c (i a +^B a . 



(20a) 
(20b) 



Equation (20i) is the 1+log slicing condition and 



Eqs. (20d,c) are the Gamma-driver shift condition. The 
extra variable B a is a contravariant vector with no den- 



sity weight. On the right-hand side of Eq. (20 :) the term 
dt A° is eliminated using the BSSN equation of motion for 



IV. LAGRANGIAN AND EULERIAN CASES 

For the Lagrangian case we have d t g — and d t A = 0. 
For any choice of initial values, g and A will remain un- 
changed throughout the evolution. Let us call these ini- 



tial values g and A. Thus g is a time independent spatial 
scalar density of weight +2, and A is a time independent 
spatial scalar. The time independence of g and A imply 
d x In g = -2D a {3 a and d±A = -f3 a d a A. 

Since g and A equal their initial values for all time, 
we can replace g with g and A with A wherever they 
appear in Eqs. (18). Note that the covariant divergence 



d t B a = f3 c D c B a + (d t A a - /3 c D c A a ) - ^B a .(20c) 9 ^ 9 



of the shift depends on the spatial metric only through 
its determinant: D a (5 a — d a {\/g(i a ) / 'y/g. Since the de- 
terminant is constant in time, we can replace D a (3 a with 
D a (3 a , where D a is the covariant derivative built from the 
initial conformal metric g ab . If we make the replacements 
A — > A, and D a (3 a — > D a (3 a everywhere, we ob- 
tain the traditional BSSN equations written in covariant 
form. 

The traditional BSSN equations are not strongly hy- 
perbolic unless the algebraic constraint A = A is contin- 
uously enforced [TU]. (This is true for any choice of 
gauge conditions, not just the standard gauge.) In prac- 
tice, the constraint A = is imposed by making the re- 
placement A ab — > A ab — g ab A/3 after every sub-timestep 
in the numerical evolution. This prevents A from devel- 
oping a non-zero value due to numerical error. 

As an alternative, we can achieve strong hyperbolicity 
by leaving the term 2ag ab A/3 in Eq. (18i) alone [TU]. If 
we do this and also choose A = 0, we find the following 
Lagrangian BSSN equations: 



d±g a b 
d±A ab 

d ± K 
dxA a 



2 - 2 
^gabD c (3 c - 2aA ab + -ag ab A , 

2 

-A ab D c (3 c - 2aA ac A c b + aA ab K 



e~ Atp [-2aD a D b <p + 4aD a (pD b ip + 4D {a aD b) (p - D a D b a + aR ab ] TF 

D C (3 C -\aK , 
6 

K 2 + aA ab A ab - e- ilp (D 2 a + 2D a aD a ip) , 



C b D b f3 a + g bc D b D c (3 a - g bc R a bcd (3 d + -AT a D c fi c + -D a D c f3 c 



-2A bc (S^d c a ~ 6aS^d c ip - aATJ.) - -ag ab d b K 



(21a) 

(21b) 
(21c) 
(21d) 

(21e) 



r 



lZ ab is given by Eq. (19 1. These equations with the stan- 



dard gauge Eqs. (20) are strongly hyperbolic without the 
need to enforce ^4 = explicitly. If the term 2ag ab A/3 is 



omitted from Eq. (21 1), the result is the traditional BSSN 
equations. The traditional BSSN equations coincide with 
the equations that are in widespread use in the numerical 
relativity community when the background connection 
vanishes, T^ c = 0, and the initial data satisfies g = 1. 

It is a bit of an overstatement to say that the La- 



grangian system (21 1, or the traditional system with 



A = enforced, is strongly hyperbolic. As shown by 
Beyer and Sarbach [TTJ , the traditional BSSN system plus 
standard gauge is strongly hyperbolic for 2a ^ e Alfi . The 
condition 2a ^ e Av is likely violated on a surface of co- 
dimension one in black hole simulations, but in practice 
this does not seem to be a problem. 

Another natural choice for the evolution of g is the 
Eulerian case, d±g = 0. Let us assume as before that 
A is time independent, d t A = 0. Let us replace A with 
its initial value A = everywhere except in the d±g ab 
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equation. Then the Eulerian BSSN equations are 



d±g a b = -2aA ab + -ag ab A , (22a) 

d x A ab = -2aA ac A c b + aA ab K + e- 4lp [~2aD a D b ip + AaD ai pD b ip + AD [a aD b) ip - D a D b a + aU ab ] TF , (22b) 
d±(fi 
d ± K 

d x A a 



aK , 

6 

~K 2 + aA nl , \ 



iv (D 2 a + 2D a aD a ip) , 



C b D b p a + g bc D b D c p a - g bc R a bcd f3 d - 2A bc (S^d c a - 6a<5^ - aAT a bc ) - -ag ab d b K 



(22c) 
(22d) 
(22e) 



Again, lZ ab is given by Eq. (19 1. These equations with 



the standard gauge Eqs. (20) are strongly hyperbolic for 



2a ^ e Atp and do not require enforcement of the algebraic 
constraint A = A = 0. They are more simple than both 
the Lagrangian equations (21) and the traditional BSSN 
equations. 



Note that the Gamma driver shift (20:>,c) depends on 
dtA a . This term is to be replaced with appropriate terms 
from the equation of motion for A a . There are two pos- 
sibilities. The term i9 4 A a can be defined using either 



the Lagrangian equation (21 3) or the Eulerian equation 
p2^) . If one wants the Gamma driver shift condition as 
it is currently defined in the numerical relativity commu- 
nity, then the Lagrangian equation should be used. This 
is the case even if one chooses to evolve the BSSN vari- 
ables using the Eulerian Eqs. (22 1. The BSSN equations, 



that uses the Lagrangian equation for c? t A a to define the 
Gamma driver shift, is strongly hyperbolic for 2a ^ e 4,p . 

It would be interesting to investigate the properties of 
the shift condition defined by using the Eulerian equa- 
tion ( 22 3) to eliminate d t A a from Eqs. (pp|. In this case 



either Lagrangian or Eulerian, with the standard gauge 



the BSSN equations (either Lagrangian or Eulerian) plus 
gauge conditions are strongly hyperbolic for 8a ^ 3e 4<p . 
A detailed analysis of hyperbolicity for these systems will 
be given elsewhere [TU] . 
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